ZHe data
Data from Pikes Peak, Mt. Evans, Longs Peak
Rachel Havranek
28 Feb 2021
AFT_data <- read.csv("data/AFT_Data.csv")
model_outcomes <- read.csv ("data/Model_outcomes.csv")Pikes Peak
Load and Manipulate Data
Pikes_Data <- read.csv("data/Pikes_Data_for_R .csv")
#Modify the sample names, for ease of reading, where PP still stands for Pikes Peak, and they are numbered 1 - 6 with the lowest elevation sample being 1 and highest elevation sample being 6
Pikes_Data <- Pikes_Data %>%
mutate (papername =
ifelse(Sample == "PP2084", "PP1", #If True, add label PP2
ifelse( Sample == "PP2479", "PP2", #If True, add label PP2
ifelse (Sample == "PP2907", "PP3", #If true, add label PP3
ifelse (Sample == "PP3597", "PP4", #If true, add label PP4
ifelse (Sample == "PP3971", "PP5", "PP6") #if true, add label PP5, ELSE add the label PP6
)
)
)
)
)
Pikes_Data <- Pikes_Data %>% mutate(Peak = "Pikes Peak")Plots
This section makes plots for:
1. Pikes Peak Elevation - date
Pikes Peak Date-eU
May 01, 2020 – Becky and I decided to use a grey scale for the lowest three elevation samples for the Pikes Peak samples. We chose to do this because those data were previously published in Flowers et al., 2020 PNAS. I have not changed the code here, but instead I have just manually changed the colors in illustrator. I chose this approach becuase it is very hard to add two different color scales to a ggplot (they’ve done that on purpose).- Pikes Peak Elevation vs. date with each point colored by eU
- Pikes Peak date vs. grain size
Pikes Peak AFT data (Kelley and Chapin, 2004)
## [1] "Excludes bottom two AFT samples at 1777m and 1866m"## Warning: Removed 2 rows containing missing values (geom_point).
Hefty Inputs
Grains excluded from models:
* Sample: PP1, grain:Z32, rownumber: 6 - this grain has an eU of 151.5 and a date of 312.9, which is ~ 300 Ma younger than grains w/ comprable eU
* Sample: PP6, grain: Z3, rownumber: 30 - this grain had an eU of 802.5 and a date of 414, which an outlier from the rest of the bin.
We have setaparated the data into two groups for modeling - below 3500 m, and above 3500 m. This separates the data that were previously published by Flowers et al., 2020 PNAS, with the new data presented here. We have chosen the same elevation cutoff for separting the three datasets for all three peaks. The bins I have here mimic those in Flowers et al., 2020 PNAS.
Bins:
* 0 - 250 ppm (this previously was 0 - 100 + 100 - 250, but the kintetic model had a hard time with a positive Date - eU trend) * 250 - 600 ppm This is different than before so that we don’t force a positive date-eU trend from the 0 - 250, 250 - 500 bins for the high elevation model * 500 - 1000 ppm
* 1000 - 2500 ppm
Previously tried data binning methods:
* 0 - 250 ppm (this previously was 0 - 100 + 100 - 250, but the kintetic model had a hard time with a positive Date - eU trend) * 250 - 500 ppm
* 500 - 1000 ppm
* 1000 - 2500 ppm
This approach was disgarded because we decided to model the data in two groups
I have divided by data into 5 bins with roughly equal number of grains: * 0 - 150 ppm (7 grains) * 150 - 350 ppm (6 grains) * 350 - 500 ppm (6 grains) * 500 - 900 ppm (7 grains) * 900 - 2500 ppm (5 grains)
I have chosen to model all of the samples together because they form a singular date-eU trend. I like the bins I have chosen here becuase I think they accurately capture the overall date - eU span. - 670.04, eU = 102.46
- 627.27, eU = 250.52
- 288.52, eU = 447.82
- 264.16, eU = 676.41
- 140.76, eU = 1323.12
grains.not.modeled <- c(6, 30)
Pikes_Data <- Pikes_Data %>%
mutate(
Rownumber= row_number(),
Donotuse = (Rownumber %in% grains.not.modeled),
bindata= cut(eU, c(0, 250,600,1000,2500)) #these are my bin cutoffs
) Here I calculate the parapmeters I will need for my hefty model input. This automatically saves to a CSV in my data_output folder, and can be easily accessed and shared
## # A tibble: 4 x 16
## bindata U Th Sm RawDate_mean Rawdate_15perce… Rawdate_SD mean_rs
## <fct> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 (0,250] 126. 76.3 2.3 536. 80.3 23.5 65.4
## 2 (250,6… 355. 290. 9.6 332. 49.9 129. 63.1
## 3 (600,1… 652. 404. 9.2 137. 20.5 88 114.
## 4 (1e+03… 1450. 731. 13.4 128. 19.3 60.7 75.9
## # … with 8 more variables: CorrDate_mean <dbl>, CorrDate_15percent <dbl>,
## # CorrDate_SD <dbl>, model_uncertainty <dbl>, N <int>, eU <dbl>, He <dbl>,
## # FT <dbl>## # A tibble: 4 x 16
## bindata U Th Sm RawDate_mean Rawdate_15perce… Rawdate_SD mean_rs
## <fct> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 (0,250] 99.1 62.1 1.3 565. 84.8 33.6 68.1
## 2 (250,6… 409. 196. 2.1 322. 48.2 75.9 75.1
## 3 (600,1… 816 416. 10.8 134. 20.1 16.2 95
## 4 (1e+03… 1059 662. 25.8 99 14.8 5.7 107.
## # … with 8 more variables: CorrDate_mean <dbl>, CorrDate_15percent <dbl>,
## # CorrDate_SD <dbl>, model_uncertainty <dbl>, N <int>, eU <dbl>, He <dbl>,
## # FT <dbl>Hefty Output plots
This is a dynamic plot that can be zoomed in on, and if you hover over a data point, you can see what grain it is
This is a dynamic plot that can be zoomed in on, and if you hover over a data point, you can see what grain it is
Mount Evans
Load and Manipulate Data
#Read in data (same format, with modified headers, as what will be published)
Mount_Evans <- read.csv("data/Evans Data.csv")
#Add a column with the sample names from original collection.
Mount_Evans <- Mount_Evans %>% mutate(
collection.name = ifelse (Sample == "ME1_2872", "ME10-16",
ifelse (Sample == "ME2_3596", "ME8-16",
ifelse (Sample == "ME3_3978", "ME3-16",
"ME1-16")
)
)
)
#Add in a column with sample elevation
Mount_Evans <- Mount_Evans %>% mutate(
Elevation = ifelse (Sample == "ME1_2872", 2872,
ifelse (Sample == "ME2_3596", 3596,
ifelse (Sample == "ME3_3978", 3978,
4345)
)
)
)
Mount_Evans <- Mount_Evans %>% mutate(Peak = "Mount Evans")Plots
## [1] "Excludes the two highest eU data points"
Hefty Inputs
Here I have separated the Mount Evans data into below 3500 m, and above 3500 m. I separted it like this becuase it looks like the bottom sample (ME1) has a different date-eU trend than the other three samples. Additonally, Bryant and Naeser (1980) documented a ‘roll-over’ to older AFT dates at an elevation of approximatley 3500 m, and so it reasonable to expect that samples at low elevation experienced a slightly different thermal history than the high elevation samples.
ME1 Bins:
*0 - 450 ppm - we opted to combine the first two bins because we were having a hard time getting any good fits, and it was the middle bin that was hard to match 450 - 1000 ppm
Previously: 0 - 250 ppm
250 - 500 ppm
*500 - 1000 ppm
ME2, ME3, ME4 Bins:
Excluded the two high-eU grains, becuase the kineteic model can’t really deal with those data points
0 - 450 ppm – we slightly altered the boundary of this bin as compared to all other bins to capture the grain from ME4, with an eU of 489.1 ppm and a date of 363.0 which pretty clearly falls into the higher eU population 450 - 1000 ppm 1000 - 2000
Older separation: Bins:
* 0 - 400 (9 grains)
* 400 - 700 (8 grains)
* 700 - 1600 (5 grains)
* 1600 - 5100 (2 grains)
Hefty outputs for the lowest sample:
## # A tibble: 2 x 16
## bindata U Th Sm RawDate_mean Rawdate_15perce… Rawdate_SD mean_rs
## <fct> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 (0,450] 174. 99.1 1 362. 54.3 44.6 65.4
## 2 (450,1… 637. 303. 2.3 136. 20.4 38.3 86.8
## # … with 8 more variables: CorrDate_mean <dbl>, CorrDate_15percent <dbl>,
## # CorrDate_SD <dbl>, model_uncertainty <dbl>, N <int>, eU <dbl>, He <dbl>,
## # FT <dbl>grains.not.modeled <- c(13, 15) #removes grains with an eU of 5086.6 ppm and 3871.2 ppm, respectively
Mount_Evans <- Mount_Evans %>%
mutate(
Rownumber= row_number(),
Donotuse = (Rownumber %in% grains.not.modeled),
bindata= cut(eU, c(0,450, 1000,2000)) #these are my bin cutoffs
) Here I calculate the parapmeters I will need for my hefty model input. This automatically saves to a CSV in my data_output folder, and can be easily accessed and shared
## # A tibble: 3 x 16
## bindata U Th Sm mean_rs RawDate_mean Rawdate_15perce… Rawdate_SD
## <fct> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 (0,450] 294. 116. 3 69.1 433. 64.9 51.6
## 2 (450,1… 617. 265. 4.8 72.5 253. 38 69.7
## 3 (1e+03… 1195. 188. 9.9 62.6 141. 21.2 36.6
## # … with 8 more variables: CorrDate_mean <dbl>, CorrDate_15percent <dbl>,
## # CorrDate_SD <dbl>, model_uncertainty <dbl>, N <int>, eU <dbl>, He <dbl>,
## # FT <dbl>Hefty outputs
This is a dynamic plot that can be zoomed in on, and if you hover over a data point, you can see what grain it is
Longs Peak
Load and Manipulate Data
Longs_Peak <- read.csv("data/Longs_Peak.csv")
#Rename sample to match other peak schema
# Longs_Peak <- Longs_Peak %>% mutate (
# papername = ifelse (Sample == "LP1", "LP7",
# ifelse (Sample == "LP2", "LP6",
# ifelse (Sample == "LP3", "LP5",
# ifelse (Sample == "LP4", "LP4",
# ifelse (Sample == "LP7", "LP3",
# ifelse(Sample == "LP5", "LP2", "LP1")
# )
# )
# )
# )
# )
# )
# Add in column with sample elevation
Longs_Peak <- Longs_Peak %>% mutate (
Elevation = ifelse (Sample == "LP7", 4343,
ifelse (Sample == "LP6", 4121,
ifelse (Sample == "LP5", 4023,
ifelse (Sample == "LP4", 3688,
ifelse (Sample == "LP3", 3500,
ifelse(Sample == "LP2", 3383, 2835)
)
)
)
)
)
)
Longs_Peak <- Longs_Peak %>% mutate (Peak = "Longs Peak")Plots
Hefty Outputs
Just like the other two peaks, I’ve divided this dataset into a low elevation (<= 3500 m), and high elevation (>3500 m) groups.
Even eU bins - I think this makes the most sense in some ways for this distribution, because it honors ‘clusters’ for the most part, and evenly divides the eU space. Downside: it is not particularly similar to either of the groupings from the other peaks 4 Bins: 0 - 500 500 - 1000
*1000 - 2200 None of these options really meaningfully change the average dates or the uncertainties (either 15% or s.d.) in each resepctive bin.
grains.not.modeled <- c()
Longs_Peak <- Longs_Peak %>%
mutate(
Rownumber= row_number(),
Donotuse = (Rownumber %in% grains.not.modeled),
bindata= cut(eU, c(0,500, 1000, 3200)) #these are my bin cutoffs
) Here I calculate the parapmeters I will need for my hefty model input. This automatically saves to a CSV in my data_output folder, and can be easily accessed and shared
## # A tibble: 3 x 15
## bindata U Th RawDate_mean Rawdate_15perce… Rawdate_SD mean_rs
## <fct> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 (0,500] 452. 116. 99.6 15.0 NA 43.3
## 2 (500,1… 708. 87.4 36.3 5.44 5.84 39.0
## 3 (1e+03… 1446. 452. 38.7 5.81 8.98 37.8
## # … with 8 more variables: CorrDate_mean <dbl>, CorrDate_15percent <dbl>,
## # CorrDate_SD <dbl>, model_uncertainty <dbl>, N <int>, eU <dbl>, He <dbl>,
## # FT <dbl>## # A tibble: 3 x 15
## bindata U Th RawDate_mean Rawdate_15perce… Rawdate_SD mean_rs
## <fct> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 (0,500] 374. 51.2 54.8 8.23 NA 50.9
## 2 (500,1… 750. 161. 46.9 7.03 8.21 45.8
## 3 (1e+03… 2442. 255. 34.0 5.1 2.02 38.0
## # … with 8 more variables: CorrDate_mean <dbl>, CorrDate_15percent <dbl>,
## # CorrDate_SD <dbl>, model_uncertainty <dbl>, N <int>, eU <dbl>, He <dbl>,
## # FT <dbl>Hefty outputs
Other Plots
CL_raman_comparison <- Pikes_Data %>%
filter(Sample == "PP2084") %>%
ggplot()+
geom_point(aes(eU, Corr_Date, color = Grain, label = Grain))## Warning: Ignoring unknown aesthetics: labelCL_raman_comparison
date_overlap_highelevation <- Pikes_Data %>% filter(papername == "PP4" | papername == "PP5" | papername == "PP6") %>%
bind_rows(
Mount_Evans %>% filter (Sample == "ME2_3596" | Sample == "ME3_3978" | Sample == "ME4_4345")
) %>%
ggplot()+
aes(x = eU, y = Corr_Date, shape = Peak)+
scale_shape_manual(values=c(15, 17)) +
geom_point(size = 3)+
labs (title = "High Elevation Samples")## Warning in bind_rows_(x, .id): Unequal factor levels: coercing to character## Warning in bind_rows_(x, .id): binding character and factor vector, coercing
## into character vector
## Warning in bind_rows_(x, .id): binding character and factor vector, coercing
## into character vector## Warning in bind_rows_(x, .id): Unequal factor levels: coercing to character## Warning in bind_rows_(x, .id): binding character and factor vector, coercing
## into character vector
## Warning in bind_rows_(x, .id): binding character and factor vector, coercing
## into character vector## Warning in bind_rows_(x, .id): Unequal factor levels: coercing to character## Warning in bind_rows_(x, .id): binding character and factor vector, coercing
## into character vector
## Warning in bind_rows_(x, .id): binding character and factor vector, coercing
## into character vectordate_overlap_highelevation
date_overlap_lowelevation <- Pikes_Data %>% filter(papername == "PP3") %>%
bind_rows(
Mount_Evans %>% filter (Sample == "ME1_2872")
) %>%
ggplot()+
aes(x = eU, y = Corr_Date, shape = Peak)+
scale_shape_manual(values=c(15, 17)) +
geom_point(size = 3)+
labs (title = "Low Elevation Samples")## Warning in bind_rows_(x, .id): Unequal factor levels: coercing to character
## Warning in bind_rows_(x, .id): binding character and factor vector, coercing
## into character vector
## Warning in bind_rows_(x, .id): binding character and factor vector, coercing
## into character vector## Warning in bind_rows_(x, .id): Unequal factor levels: coercing to character## Warning in bind_rows_(x, .id): binding character and factor vector, coercing
## into character vector
## Warning in bind_rows_(x, .id): binding character and factor vector, coercing
## into character vector## Warning in bind_rows_(x, .id): Unequal factor levels: coercing to character## Warning in bind_rows_(x, .id): binding character and factor vector, coercing
## into character vector
## Warning in bind_rows_(x, .id): binding character and factor vector, coercing
## into character vectordate_overlap_lowelevation